NAVAL  POSTGRADUATE  SCHOOL 
Monterey,  California 


THESIS 


Approved  for  public  release;  distribution  is  unlimited. 


20020102  091 


REPORT  DOCUMENTATION  PAGE  |  Form  Approved  OMB  No.  0704-0188 

Public  reporting  burden  for  this  collection  of  information  is  estimated  to  average  1  hour  per  response,  including 
the  time  for  reviewing  instruction,  searching  existing  data  sources,  gathering  and  maintaining  the  data  needed,  and 
completing  and  reviewing  the  collection  of  information.  Send  comments  regarding  this  burden  estimate  or  any 
other  aspect  of  this  collection  of  information,  including  suggestions  for  reducing  this  burden,  to  Washington 
headquarters  Services,  Directorate  for  Information  Operations  and  Reports,  1215  Jefferson  Davis  Highway,  Suite 
1204,  Arlington,  VA  22202-4302,  and  to  the  Office  of  Management  and  Budget,  Paperwork  Reduction  Project 
(0704-0188)  Washington  DC  20503.  _ 


1.  AGENCY  USE  ONLY  (Leave  blank)  1  2.  REPORT  DATE  1  3.  REPORT  TYPE  AND  DATES  COVERED 

I  June  2001  I  Master’s  Thesis 


.  TITLE  AND  SUBTITLE:  AUV  Steering  Parameter  Identification  For  5.  FUNDING  NUMBERS 
roved  Control  Design  M000 1 400AF00002 

6.  AUTHOR(S)  Jay  H.  Johnson _ _ 


7.  PERFORMING  ORGANIZATION  NAME(S)  AND  ADDRESS(ES) 

Naval  Postgraduate  School 

Monterey,  CA  93943-5000 _ 


9.  SPONSORING  /  MONITORING  AGENCY  NAME(S)  AND  ADDRESS(ES) 
Office  of  Naval  Research,  800  N.  Quincy  St.,  Arlington,  VA  22217- 
5660 


10.  SPONSORING  /  MONITORING 
AGENCY  REPORT  NUMBER 


11.  SUPPLEMENTARY  NOTES  The  views  expressed  in  this  thesis  are  those  of  the  author  and  do  not  reflect  the  official 
policy  or  position  of  the  Department  of  Defense  or  the  U.S.  Government. _ 


12a.  DISTRIBUTION  /  AVAILABILITY  STATEMENT  I  12b.  DISTRIBUTION  CODE 

Approved  for  public  release:  distribution  is  unlimited.  I 


13.  ABSTRACT 

Any  effort  to  provide  precision  control  for  an  Autonomous  Underwater  Vehicle  requires  an  accurate 
estimation  of  both  the  vehicles  physical  and  hydrodynamic  parameters.  Here  a  vehicle  model  for  controlled 
steering  behaviors  was  developed  and  the  hydrodynamic  parameters  were  calculated  from  actual  data  obtained 
from  operation.  The  steering  equation  parameters  are  based  on  a  least  squares  fit  to  sideslip  and  turn  rate  data 
using  maximum  likelihood  of  batch  processing.  In  this  way,  a  more  accurate  simulation  has  been  found  for  the 
development  of  a  track  controller  that  stably  drives  the  vehicle  between  mission  waypoints.  Prediction  accuracy 
of  the  model  was  better  than  ninety-five  percent  over  the  data  set  used. 


14.  SUBJECT  TERMS  Underwater  Vehicle,  AUV,  Control,  System  Identification 


17.  SECURITY 
CLASSIFICATION  OF 
REPORT 

Unclassified 


18.  SECURITY 
CLASSIFICATION  OF  THIS 
PAGE 

Unclassified 


19.  SECURITY 
CLASSIFICATION  OF 
ABSTRACT 

Unclassified 


15.  NUMBER  OF 
PAGES  ^ 

72 


16.  PRICE  CODE 


20.  LIMITATION 
OF  ABSTRACT 


THIS  PAGE  INTENTIONALLY  LEFT  BLANK 


11 


Approved  for  public  release;  distribution  is  unlimited. 


AUV  STEERING  PARAMETER  IDENTIFICATION  FOR  IMPROVED 

CONTROL  DESIGN 

Jay  H.  Johnson 

Lieutenant,  United  States  Navy 
B.S.,  University  of  Kansas,  1993 


Submitted  in  partial  fulfillment  of  the 
requirements  for  the  degree  of 


MASTER  OF  SCIENCE  IN  MECHANICAL  ENGINEERING 


from  the 


NAVAL  POSTGRADUATE  SCHOOL 
June  2001 


Author:  - 

'  ^  Jay  H.  Johnson 


111 


THIS  PAGE  INTENTIONALLY  LEFT  BLANK 


IV 


ABSTRACT 


Any  effort  to  provide  precision  control  for  an  Autonomous  Underwater  Vehicle 
requires  an  accurate  estimation  of  both  the  vehicles  physical  and  hydrodynamic 
parameters.  Here  a  vehicle  model  for  controlled  steering  behaviors  was  developed  and 
the  hydrodynamic  parameters  were  calculated  from  actual  data  obtained  from  operation. 
The  steering  equation  parameters  are  based  on  a  least  squares  fit  to  sideslip  and  turn  rate 
data  using  maximum  likelihood  of  batch  processing.  In  this  way,  a  more  accurate 
simulation  has  been  found  for  the  development  of  a  track  controller  that  stably  drives  the 
vehicle  between  mission  waypoints.  Prediction  accuracy  of  the  model  was  better  than 
ninety-five  percent  over  the  data  set  used. 
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I.  INTRODUCTION 


A.  GENERAL  BACKGROUND  AND  LITERATURE 

The  increased  desire  to  use  Autonomous  Underwater  Vehicles  (AUVs)  and 
Unmanned  Undersea  Vehicles  (UUVs)  for  commercial  and  military  applications  has  led 
to  a  great  deal  of  research  in  this  field  over  the  last  decade.  Specifically,  the  military 
need  for  clandestine  mine  reconnaissance  is  now  being  solved  using  AUVs.  A  high 
degree  of  maneuverability  in  shallow  water  conditions  means  these  vehicles  must  have 
high  performance  steering  controllers  requiring  accurate  computer  models.  Although 
computer  models  have  facilitated  the  estimation  of  fluid  parameters  on  vehicles,  the 
building  of  accurate  models  to  match  complex  shapes  is  tedious  and  time  consuming.  An 
exact  calculation  of  the  dynamic  parameters  of  an  AUV  is  impossible.  Even  with  an 
estimation  of  these  parameters,  there  remains  the  question  of  validation  of  the  results 
obtained.  To  get  the  hydrodynamic  properties  for  marine  vehicles  generally  one  of  the 
three  following  methods  is  normally  used: 

1 .  Detailed  computer  modeling  as  previously  mentioned; 

2.  A  tow  tank.  In  a  tow  tank,  generally  a  model  (models  of  large  vessels  are 
used  due  to  size  restrictions)  of  the  vehicle  is  pulled  through  the  water  and 
maneuvered  while  maintaining  a  constant  speed.  During  this  entire  time, 
instruments  are  used  to  measure  the  forces  acting  on  the  vehicle.  This  process  is 
very  informative  and  provides  very  good  results.  But,  it  is  very  time  consuming 
and  expensive  to  accomplish;  and 

3.  By  approximating  the  shape  of  the  vehicle  to  known  shapes  and  using  this 
approximation  to  calculate  some  coefficients  and  using  system  identification 
methods  for  others. 
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Apart  from  the  calculation  of  the  hydrodynamic  masses,  hydrodynamic 
coefficients,  Yv,  Yr,  Nv,  Nr,  Yg  and  Ng  are  not  easily  computed,  and  when 
measured  using  2,  do  not  always  match  field  data  from  3. 

B.  SCOPE  OF  THIS  WORK 

The  overall  problem  of  marine  vehicle  hydrodynamics  is  very  complex  and 
diverse.  This  study  will  focus  on  the  identification  of  the  NPS  ARIES  AUV  steering 
system  hydrodynamic  parameters.  In  particular  the  values  of  added  mass  for  sway  and 
yaw  are  determined  from  physical  measurements,  which  are  adapted  to  an  approximate 
shape.  Normal  system  identification  techniques  [Ljung  et  al.]  obtain  values  for  a 
canonical  set  of  transfer  function  parameters  for  which  Tv,  Yr,  Nv  and  Nr  (defined  in 
Chapter  II)  cannot  be  extracted  uniquely.  However,  with  turn  rate  (r),  sideslip  velocity 
(v),  and  rudder  angle  (S)  measurements,  it  is  possible  to  get  a  uniquely  identifiable 
estimate  of  each  parameter  Nv,  Nr,  Yv  and  Yr.  In  this  case,  we  use  a  set  of  data  obtained 
from  ARIES  and  use  a  least  mean  square  batch  regression  process  to  obtain  unique 
values. 

Previous  literature  on  this  subject  includes  a  major  work  by  Abkowitz,  (1980) 
using  a  tanker.  In  which  the  estimated  hydrodynamic  coefficients  were  obtained  using  an 
extended  Kalman  Filter.  The  results  were  inconclusive  however,  because  they  did  not 
match  model  test  data  and  the  ship  had  to  be  built  in  order  to  get  test  results. 

Chapter  II  discusses  the  equations  of  motion  for  an  underwater  vehicle  and  AUV 
modeling.  It  also  includes  the  calculations,  formulas  and  approximations  used  and  the 
reasoning  behind  each. 
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Chapter  III  discusses  the  theory  of  parameter  identification  from  vehicle 
measurements  and  the  methodology  used  in  this  report. 

Chapter  IV  provides  information  on  the  data  collected  by  the  ARIES  AUV  and 
the  processing  required  to  obtain  the  final  results. 

Chapter  V  proves  that  the  methodology  used  in  this  report  provides  an  accurate 
way  to  obtain  hydrodynamic  coefficients  from  vehicle  data. 

Chapter  VI  lists  the  conclusions  of  this  report  derived  from  Chapter  V. 
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II.  EQUATIONS  OF  MOTION  AND  AUV  MODELING 


A.  GENERALIZED  STEERING  EQUATIONS  OF  MOTION 

This  section  will  discuss  the  equations  of  motion  for  a  rigid  body  in  an  underwater 
environment  leading  to  the  model  used  for  AUV  steering  control  design  and  the 
simplifications  made  for  the  analysis  of  the  NPS  AUV. 

The  forces  acting  on  the  AUV  determine  its  motion  and  are  independent  of  the 
inertial  position.  The  equations  of  motion  are  given  using  Euler  equations  for  the  body 
fixed  frame  (Greenwood,  1980).  The  six  simplified  equations  of  motion  given  by 
(Healey,  1995)  are: 

SURGE  EQUATION  OF  MOTION 

m[ur-vrr  +  wrq-xG(q 2  +r2)+  yG(pq-r)+za(pr  +  q)]  +(iV- #)sin  0  =  Xf  (1) 

SWAY  EQUATION  OF  MOTION 

m[  vr  +  urr  -  wr  p  +  xG  (pq + r)-yG  (p2  +  r 2 )+  zG  (qr  -  />)]  -(w  -  B)  cos  0  sin  <|>  =  Yf  (2) 

HEAVE  EQUATION  OF  MOTION 

m[wr-urq  +  vrp  +  xG(pr-q)+yG(qr  +  p)-zG[p2 +q2)]+{W-B)cosQcos$  =  Z/  (3) 

ROLL  EQUATION  OF  MOTION 

lx P  +  (a  -  ■ I y  V  +  ^ hy  ipr  -  i)~  I yz  (v‘ 1  “  r  2  )“  4z  (w  +  ^ ^)+  G  ~ « r 9' + vr  p)  (4) 

—  zG(vr  +urr -wrp)]-(yGW  -j'jSJB)cos0cos<|>  +  (zGfF'  -z55)cos0sini|)  =  K  j 
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PITCH  EQUATION  OF  MOTION 


^  +  (A  “ A )pr  ~!xy .(qr  +  phlyzipq-rj+Ixzip2 -r2)-m[xc(w-urq+vrp)  (5) 

—  Zq (ur  -vrr+wrq)]+(xGW -xbB)cosQcos$  +  (zgW -zBB)sinQ  =  M j 

YAW  EQUATION  OF  MOTION 

I2r  +  (ly  -  Ix )pq  -I^p2  -q2)- I>z  (pr  +  q)+  (qr  -p)+  m[xG (vr  +  urr  -  wr p)  (6) 

-ycfar  ~  vrr  +  wr?)]-(xG  W-xB  B)  cos  ©sin  <J)  —  (yG  W  -  yBB)s  in  6  =  Nf 

The  variables  u,  v  and  w  are  the  component  velocities  expressed  in  terms  of  a 
body  fixed  coordinate  system  and  the  r  subscript  denotes  that  the  velocity  is  with  respect 
to  the  water.  The  variables  p,  q  and  r  are  the  components  of  angular  velocity  expressed 
in  the  body  fixed  coordinate  system.  W  and  B  are  weight  and  buoyancy  and  I  is  a 
moment  of  inertia  term.  The  x.  y  and  z  variables  with  the  B  and  G  subscript  come  from 
the  position  difference  between  the  center  of  gravity,  center  of  buoyancy  and  the  center 
of  the  vehicle  and  the  effects  of  gravity  and  buoyancy.  All  of  the  right  hand  terms  are  the 
forces  acting  on  the  AUV  acting  in  the  associated  direction. 

Further  simplifications  can  be  made  to  Equations  1  thru  6  with  the  following 
assumptions  for  the  AUV : 

1 .  That  the  vehicle  is  neutrally  buoyant,  B  =  G. 

2.  Horizontal  plane  motion,  [w,  ,p,q,Z,<f>,OpcG,yc\  =  0. 

3.  The  forward  speed  ur  is  a  constant  U0. 
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4.  Surge  motion  changes  are  neglected. 


5.  Roll  coupling  is  neglected. 

These  assumptions  reduce  the  steering  equations  of  motion  to: 


ur  =U0 

(7) 

mvr  =-mU0r  +  AYf{t) 

(8) 

% 

il 

a 

(9) 

\j/  =  r 

(10) 

X  -  U0  cos\|/-vr  sinvj z+t/^ 

HD 

Y -U 0  sin\j/-vr  cos\|/  +  C/^ 

(12) 

Now  that  the  equations  of  motion  have  been  reduced,  the  forces  acting  on  the 
AUV  must  be  addressed.  Recognizing  that  under  steady  motion  there  is  a  balance 
between  lift,  drag,  buoyancy,  weight,  and  propulsion  provides  the  starting  point  to 
develop  the  components  AYf(t)  and  ANf(t).  The  variables  AYf(t)  and  an,  hi  represent 
forces  as  a  function  of  the  vehicles  dynamic  parameters.  Equation  13  assumes  that  all 
motion  variables  to  be  small  compared  to  vehicle  speed  and  the  J implies  a  ‘small  change 
in’. 

5F  =  f( ur,  vr,  wr,p,  q,  r,  du/dt,  dv/dt,  w/dt,  dp/ dp,  dq/dt,  dr/dt,  t)  (13) 

Simplifying  Equation  13  is  based  on  the  following  rational  that  the  AUV  is 
generally  symmetric  about  its  longitudinal  axis  and  that  fluid  flow  in  one  direction  will 
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not  greatly  affect  an  orthogonal  direction.  The  expression  for  the  transverse  force 
becomes: 


Yf=Ytvr+Yvvr+Yff  +  Yrr 


(14) 


and 


n  = 


dIf_-  Y  •  v  _£Xl- Y  -aYf- 

a..  5  1 V.  —  o.  5  1  r  —  5  1  r  ~ " 

ovr  5vr  5r 


dr  ’ 


evaluated  at  the  forward  speed  U0.  The  coefficient  is  the  added  mass  in  sway, 
evaluated  with  the  vehicle  in  forward  motion,  Yi  is  added  mass  induced  by  yaw,  YVr  is  a 
coefficient  of  sway  force  from  side  slip,  and  Yr  is  a  coefficient  of  sway  force  from  yaw. 
The  rudder  when  actuated  will  produce  a  lift  force  that  when  linearized  is  expressed  as 
Ys5r(t).  Likewise  the  rotational  force  is  developed  as: 

Nf  =NirVr+NVrVr+Ni.r  +  Nrr  (15) 


and 


=  ^_  .at.  N.  JNt.. 

r  5vr  r  dvT  dr 


Nr  = 


SNf  . 

~aT’ 


evaluated  at  the  forward  speed  U0.  The  coefficient  Nir  is  the  added  mass  moment  of 
inertia  in  sway,  evaluated  with  the  vehicle  in  forward  motion,  Nfis  added  mass  moment 
induced  by  yaw,  NVr  is  a  coefficient  of  yaw  moment  from  side  slipping,  and  Nr  is  a 

coefficient  of  sway  moment  from  yaw.  The  rudder  when  actuated  will  produce  a  force 
that  linearized,  is  N58r(t) .  The  equations  of  motion  in  the  horizontal  plane  are  given  by: 


Izz  r  =  Nir  vr  +  NVr  v  r  +  Nj.r  +  Nrr  +Nhbr  (t) 


(17) 


When  the  equations  are  cast  into  matrix  form  of  Mx  =Ax  +  Bn  the  result  is 


YVr  Yr  -mU0~ 
NVr  Nr 


(18) 


In  Equation  18  the  term  8,(/)  is  a  generalized  command  that  represents  the  control 
input  to  both  rudders.  It  should  be  noted  that  while  the  rudders  act  together  they  turn  in 
opposite  directions.  This  tandem  but  opposite  operation  allows  for  rapid  turning  during 
operation. 

The  final  approximation  made  for  the  AUV  was  that  the  cross  coupling  terms  in 
the  mass  matrix  are  zero.  This  approximation  is  based  on  the  fact  that  the  vehicle  is 
symmetric  and  the  rudders  are  very  close  to  being  equidistant  from  the  body  center. 
Therefore  the  final  set  of  equations  to  be  solved  for  steering  behaviors  are: 


The  model  is  complete  if  values  for  the  hydrodynamic  coefficients  can  be  found, 
and  is  appropriate  for  low  angle  of  sideslip  maneuvering  in  a  longitudinal  plane. 

B.  AUV  GEOMETRY,  INERTIA  PARAMETERS 

The  usefulness  of  any  model  simulation  is  directly  linked  to  the  accuracy  of  the 
model  being  used.  In  theory  an  exact  model  will  give  an  exact  prediction  of  results  for  a 
particular  simulation.  However,  in  reality  a  perfect  model,  except  in  limited  cases  can 
never  be  achieved.  The  purpose  of  this  section  is  to  describe  the  NPS  ARIES  AUV,  how 
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it  was  modeled  and  discusses  assumptions  and  approximations  that  were  made  to  develop 
the  final  model  used. 

1.  AUV  Model  General  Parameters 

The  AIRES  AUV  shown  in  Figure  1  is  the  vehicle  being  modeled  in  this  work. 
The  parameters  required  for  analysis  include  weight,  volume,  moment  of  inertia  and 
hydrodynamic  added  mass.  To  begin  modeling  the  AUV  the  overall  dimensions  were 
measured.  Measured  dimensions  obtained  were  as  follows: 

1 .  Length  from  nose  to  end  of  tail  section  (not  including  thrusters)  1 2$  inches. 

2.  Width  determined  longitudinally  at  the  center  section  16  inches. 

3.  Height  determined  vertically  at  center  section  10  inches. 

4.  Weight  was  recorded  on  a  lifting  scale  to  be  490  lbs. 


The  vehicle  was  then  divided  into  three  distinct  sections  for  analysis.  The  tail 
section  was  modeled  as  a  triangular  prism.  The  center  section  was  modeled  as  a 
rectangular  prism.  And  finally,  the  nose  section  modeled  as  a  half  ellipsoid. 
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The  assumption  of  symmetry  for  the  AUV  defines  the  body  axes  as  the  principal 
axis  set  for  the  vehicle.  When  evaluating  the  inertia  matrix  it  then  becomes  diagonal  (I**, 

ivy-  I zz)- 

2.  AUV  Tail  Section  Analysis 


The  AIRES  AUV  tail  section  was  modeled  as  a  rectangular  prism  shown  in 
Figure  1.  The  dimensions  used  were  based  on  where  the  vehicle  started  to  change  shape 
at  the  stem.  For  this  and  all  subsequent  discussions  of  parameters  the  x-axis  is  defined  as 
the  longitudinal  or  lengthwise  direction  (length  nose  to  tail).  The  y-axis  is  defined  to  be 
the  athwartships  or  beam  direction  (width).  Lastly,  the  z-axis  represents  the  overall 
height  of  the  vehicle  in  the  vertical  direction  (height).  The  values  measured  for 
calculations  involving  the  tail  section  are: 


a.  Height 

10  inches 

b.  Width 

16  inches 

c.  Length 

26  inches 

The  calculation  of  volume  was  conducted  by  multiplying  the  area  of  the  side  face 
(-|  Length* Height)  by  the  width.  Hence  the  total  volume  of  the  tail  section  was 
determined  to  be  2080  in3. 


II 


Figure  2  Tail  Section  Model 


The  next  step  was  to  calculate  the  mass  moment  of  inertia  for  the  tail  section. 
This  required  two  steps.  The  first  was  to  determine  the  moment  of  inertia  for  the  section 
itself  and  the  second  was  to  translate  the  calculated  moment  to  a  selected  reference  point 
on  the  AUV.  All  moments  of  inertia  were  translated  to  the  geometric  center  of  the  mid 
section.  The  moment  of  inertia  for  the  tail  section  was  not  tabulated  in  Beer  and  Johnson 
(1997),  Rothbart  (1964)  or  any  of  numerous  texts  reviewed.  The  calculation  was  carried 
out  using  principals  contained  in  Beer  and  Johnson  (1997),  these  can  also  be  found  in 
most  calculus  textbooks,  i.e.  Finney  and  Thomas  (1994).  The  solution  for  the  moment  of 
inertia  in  the  x-direction  was  accomplished  as  follows: 

First  a  differential  mass  (dm)  was  chosen  using  a  strip  in  the  x-y  plane.  The 
origin  was  chosen  to  be  on  the  forward  face  center.  The  differential  equation  of  mass 


was  developed  based  on  geometry.  Performing  the  integration  of  both  sides  yields  the 
result. 


m-p—abc 

2 


or  m  =  pV 


(21) 


where  m  is  the  mass,  p  is  the  density,  a  b  and  c  are  the  dimensions,  and  V  is  the  volume. 
For  the  elemental  rectangular  area  in  the  y-z  plane,  lMrea=  • b2dm .  Thus  the 

differential  mass  moment  of  inertia  is  ptl  zarea  where  t  is  the  thickness  of  the 
differential  mass. 

Then  the  total  differential  mass  moment  of  inertia  (dlz )  is, 

dl=  =  ptl. area  Jrx2dm.  (22) 

Where  the  jrdm  term  translates  the  moment  to  the  center  of  the  AUV  mid  section 
using  the  parallel  axis  theorem  as  outlined  in  Beer  and  Johnson  (1997).  Integrating  dl~ 

from  x  =  0  to  26  inches  for  ptl  .area  and  x  =  37.5  to  63.5  inches  for  the  r?dm  term  and 
substituting  M  for  pV, 

Izz  =  2190.25*M. 

Similar  calculations  were  preformed  for  Ixx  and  I>y  using  MAPLE  software  and 
are  included  in  Appendix  A.  The  resulting  values  are  shown  in  Table  1 . 
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Table  1  Mass  Moment  of  Inertia  for  AUV  Tail  Section 


Axis 

Mass  Moment  of  Inertia 

Units 

Ixx 

142.33*M 

lb*in*sec2 

Ivy 

2177.25*M 

lb*in*sec2 

i 

Lz 

2190.25*M 

lb*in*sec2  j 

1 

3.  AUV  Mid  Section  Analysis 

The  center  section  was  treated  as  rectangular  prism  shown  in  Figure  3  whose 
dimensions  are: 

a.  Height  10  inches 

b.  Width  16  inches 

c.  Length  75  inches 


Figure  3  Mid  Section  Model 

The  formulas  describing  the  mass  moment  of  inertia  for  a  rectangular  prism  as 
shown  in  Equations  23  thru  25  are  readily  available  in  many  textbooks. 
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(23) 


I„=±M{b2+c2) 

I>y=±M(a2+c2)  (24) 

Izz  =  ±M(a2+b2)  (25) 


including  Beer  and  Johnson  (1997)  and  Finney  and  Thomas  (1994).  The  derivation  of 
the  formulas  is  omitted  and  only  the  results  are  shown  in  Table  2  below. 


able  2  Mass  Moment  of  Inertia  for  Mid  Section 


Axis 

Mass  Moment  of  Inertia 

Units 

Ixx 

29.67*M 

lb*  in* sec2 

lyy 

490.08*M 

lb*in*sec2 

hz 

477.08*M 

lb*in*sec2 

4.  AUV  Nose  Section  Analysis 

The  moment  of  inertia  calculation  was  carried  out  using  the  same  methodology  as 
the  tail  section.  The  model  used  for  the  nose  shown  in  Figure  4  is  that  of  a  half  ellipsoid. 
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Figure  4  AUV  Nose  Section  Model 


First  a  differential  mass  (dm)  was  chosen  using  a  strip  in  the  y-z  plane.  The 
origin  was  chosen  to  be  on  the  aft  face  center.  The  differential  equation  of  mass 


/ 

dm  =  p'zbc 

V 


(26) 


was  developed  based  on  geometry  of  an  ellipsoid.  Performing  the  integration  of  both 
sides  yields  the  result 


m  =  pTtabc  or  m  =  pV 

where  m  is  the  mass,  p  is  the  density,  a  b  and  c  are  the  radial  dimensions,  and  V  is  the 
volume.  For  an  elemental  elliptical  area  in  the  y-z  plane,  1:  area  -  —zy* ,  Beer  and 
Johnson  (1997).  Thus  the  differential  mass  moment  of  inertia 


dIB 


mass  =  ptl  .area 


\4  J 


(27) 
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where  dx  is  the  thickness  of  the  differential  mass.  Then  the  total  differential  mass 
moment  of  inertia  (dlz)  is, 

dl^  =  ptl2 area + x2 dm  .  (28) 

Where  the  x2dm  term  translates  the  moment  of  the  nose  to  the  center  of  the  AUV 
mid  section  using  the  parallel  axis  theorem  as  outlined  in  Beer  and  Johnson  (1997). 

Integrating  dl-z  from  x  =  0  to  27  inches  for  ptl  zarea  and  x  =  37.5  to  63.5  inches  for  the 
x2 dm  term  and  substituting  M  for  pF, 

Izz  =  2316.43*M 

Similar  calculations  were  preformed  for  Ixx  and  Iyy  using  MAPLE  software  and 
are  included  in  Appendix  A.  The  resulting  values  are  shown  in  Table  3  below. 

Table  3  Mass  Moment  of  Inertia  for  AUV  Nose  Section 
Axis  Mass  Moment  of  Inertia  Units 

Ixx  158.60*M  lb*in*sec2 

Iw  1789.85*M  lb*in*sec2 

la  2316.43*M  lb*  in*  sec2 

The  value  of  density  chosen  for  the  AUV  was  based  on  its  weight  and 
displacement.  To  determine  the  displacement,  a  calculation  was  preformed  to  determine 
the  volume  of  seawater  for  a  weight  equivalent  to  the  AUV  minus  five  pounds  for 
buoyancy.  Since  the  nose  section  had  a  free  flood  area  the  actual  displacement  of  this 
section  had  to  be  estimated.  To  adjust  the  volume  of  the  nose  section  a  percent  fill  factor 
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was  used  to  adjust  the  volume  of  the  nose  from  a  full  volume  to  a  displaced  volume. 
Using  a  fill  factor  of  forty  percent  resulted  in  a  difference  between  calculated  and 
estimated  values  of  approximately  one  percent.  A  uniform  density  was  then  calculated 
by  dividing  the  AUV  weight  by  the  calculated  displaced  volume.  All  calculations  are 
listed  in  the  MATLAB  file  auvcalc.m  in  Appendix  B. 


Table  4  shows  the  total  mass  moments  of  inertia  as  a  result  of  summing  the  v,alue 
of  all  three  sections. 

Table  4  Total  Mass  Moment  of  Inertia  for  the  ARIES  AUV 


Axis 

Mass  Moment  of  Inertia 

Units 

Ixx 

330.60*M 

lb*  in*  sec2 

Iyy 

4457.18*M 

lb*in*sec2 

Izz 

4983.76*M 

lb*  in*  sec2 

C.  AUV  ADDED  MASS  PARAMETERS 

This  section  will  discuss  the  method  used  to  obtain  the  added  mass  coefficients 
Nr  and  Yir  of  Equation  (19).  All  of  the  actual  calculations  are  contained  in  the  MATLAB 

file  auvcalc.m  in  Appendix  B.  Added  mass  is  a  concept  that  represents  the  inertial 
reaction  on  a  moving  body  in  a  dense  medium  arising  from  the  time  rate  of  change  of 
fluid  momentum.  While  many  textbooks,  such  as  White,  1999  have  tables  for  some 
general  shapes  they  do  not  have  tables  for  combinations  of  these  shapes  or  complex 
shapes.  Due  to  the  complexity  of  the  shape  of  the  AUV,  an  approximation  of  the  overall 
shape  had  to  be  made.  Initially,  an  oblate  ellipsoid  was  chosen  since  the  ability  to  specify 
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all  three  axis  dimensions  was  desired.  The  initial  reference  used  was  Lamb,  1945.  While 
the  general  description  of  how  to  solve  this  type  of  problem  was  addressed,  no  closed 
form  solution  existed.  A  report  by  Tuckerman,  1926  discussed  the  inertia  factors  of 
ellipsoidal  shaped  airships.  This  report  discusses  the  complexity  of  the  calculations  and 
the  use  of  elliptic  integrals  to  obtain  a  solution.  It  was  then  decided  to  make  a  more 
general  simplification  to  the  model.  The  ARIES  AUV  was  modeled  as  a  prolate 
ellipsoid.  This  simplification  allowed  the  use  of  tables  and  closed  form  solutions  found 
in  Lamb,  1945  to  calculate  the  added  mass  terms  needed. 

The  first  step  in  calculating  the  added  mass  coefficients  was  to  determine  the 
neutrally  buoyant  mass.  Total  shape  volume  was  calculated  from  the  measured 
dimensions.  The  total  weight  was  then  equated  to  a  volume  of  seawater.  The  purpose  of 
this  was  to  determine  the  difference  of  volumes  since  the  nose  is  a  free  flood  area. 
Results  showed  that  a  fill  factor  of  forty  percent  was  required  to  equate  the  two  volumes. 
Since  added  mass  calculations  are  based  on  displaced  volume,  the  total  volume  was  then 
equated  to  the  prolate  ellipsoid  model.  Requiring  the  volume  to  remain  the  same  for  the 
model,  it  was  decided  to  maintain  the  overall  length  and  let  the  equivalent  diameter  be 
calculated  from  the  total  volume.  The  calculated  diameter  along  with  the  length  was  used 
in  the  eccentricity  calculations.  Formulas  used  also  required  the  neutrally  buoyant  mass. 
The  AUV  was  measured  to  have  approximately  five  pounds  of  buoyancy.  Therefore,  five 
pounds  was  added  to  the  measured  weight  to  make  the  AUV  neutrally  buoyant  for 
calculations. 

The  last  step  in  calculating  the  added  mass  terms  needed  was  to  determine  I&.  The 
calculation  of  Izz  was  used  in  the  formulas  developed  in  the  previous  sections  but  with 
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minor  modifications.  Since  the  AUV  has  lateral  and  vertical  thruster  tubes  (that  were 
empty  when  data  was  obtained)  and  no  mass  moment  of  inertia,  the  effect  of  these  tubes 
was  subtracted  based  on  axis  of  rotation  for  modeling  purposes.  Table  5  lists  the  mass 
moment  of  inertia  values  calculated  using  the  MATLAB  file  auvcalc.m  listed  in  Appendix 
B. 


Table  5  Final  Mass  Moments  of  Inertia  for  the  ARIES  AUV 


Section 

Ixx 

Iyy 

Izz 

Units 

Nose 

6.74 

76.07 

98.45 

lb*  seedin'1 

Mid 

33.44 

552.25 

537.59 

lb*sec2*in'1 

Tail 

27.81 

425.44 

427.93 

lb*sec2*in'1 

Total  (Eng) 

67.99 

1053.80 

1064.00 

lb*sec2*in'1 

Total  (SI) 

7.68 

119.10 

120.26 

kg*m2 

As  a  check  of  values  calculated  for  the  mass  moment  of  inertia,  the  radii  of 
gyration  was  calculated  for  all  three  axes  and  found  to  lie  within  the  vehicle  body  at 
reasonable  points.  The  calculations  are  shown  in  the  MATLAB  file  auvcalc.m  listed  in 
Appendix  B  and  the  results  are  listed  in  Table  6  below. 


Table  6  Radii  of  Gyration  for  The  ARIES  AUV 


X  Axis 

Y  Axis 

Z  Axis 

Units 

Distance  (Eng) 

7.32 

28.81 

28.95 

in 

Distance  (SI) 

18.60 

73.20 

73.56 

cm 
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III.  PARAMETER  IDENTIFICATION  FROM  VEHICLE 

MEASUREMENTS 


A.  INTRODUCTION  AND  BACKGROUND 

The  uncertainty  of  forces  affecting  any  underwater  vehicle  motion  requires 
knowledge  of  the  methodology  of  parameter  identification.  Basically,  the  problem  is  to 
identify  the  desired  vehicle  hydrodynamic  parameters  from  a  system  model  and 
measurements  of  other  vehicle  parameters.  This  topic  is  discussed  in  control  technology 
Soderstrum  and  Ljung  (1983),  Kalman  (1969),  Gleb  (1974)  and  has  been  developed 
especially  for  dealing  with  this  problem. 

Essential  to  the  parameter  identification  is  an  accurate  model  of  the  AUV  control 
(input)  and  response  (output).  By  manipulation  of  the  equations  of  motion,  we  can  let  the 
system  parameters  become  unknowns  and  the  variables,  as  measured,  to  be  known. 


B.  GENERAL  THEORY 

To  model  the  AUV  steering  control,  a  model  was  selected  with  the  input  output 
response  of  the  form 

y(t)=  H(t)fi(t)+v(t),  (29) 

where  y(t)  is  the  primary  output,  H(t)  is  a  matrix  of  measurements  relating  the  input  to 
output,  and  0(0  being  a  parameter  vector  which  is  ideally  constant.  The  term  v(0  is  a 
zero  mean  gaussian  white  noise. 

Selecting  a  first  order  system  as  an  example,  and  temporarily  neglecting  noise  (avoiding 
differentiation  problems)  the  model  is: 

i(r)=X(w(f)-;c(f))  with 
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Substituting  and  rewriting,  the  model  becomes: 


where  y(f) 
P  and  X. 


y{t)+Xy(t)=cXu(t)  ,  or  in  matrix  form 


>(0=[-v(f)  «(')]* 


.  ^(d=[->,(d  «(0]  mi  e(d= 


is  the  derivative  of  the  measured  out  put,  P  =  cX  and  the  independent  parameters  are 


The  AUV  measurement  data  obtained  from  sensors  is  sampled  at  eight  Hertz  and 
is  therefore  not  continuous.  Since  the  data  is  sampled  at  a  discrete  interval,  a  discrete 
time  model  must  be  used.  The  discrete  time  transfer  function  of  the  above  model  is: 


cW,zW_  £*"' 


1  -  Y  a  ,z~j 


£  uj* 

7=1 


and  in  time  difference  form  the  equation  becomes: 


y,-a]yl.i-a2yl-2--  =  b1ul_]+b2ul_2+...  or  V/ 1  »-'0) 


so  that  H(t)=  ^  0(r)=  [a^,a2,...,bx,b2,..]  ■  The  first  equation 

represents  a  matrix  containing  the  system  inputs  and  outputs  and  the  second  is  again  a 
parameter  vector. 
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C.  SOLUTIONS  USING  THE  METHOD  OF  LEAST  SQUARES 

The  purpose  of  using  a  least  squares  method  is  to  estimate  the  parameters  0(r) , 

which  are  never  exactly  known,  so  that  the  difference  between  the  actual  parameter  and 
its  estimate  are  minimized.  This  difference  is  called  equation  error  and  is  given  by, 

e(t)=(y(t)-H(t^(t))  (31) 

where  §(/)  is  the  estimate  of  0(f) .  For  the  AUV  this  involves  using  the  data  obtained 
from  the  vehicles  sensors  combined  with  the  system  model  shown  in  Equation  (19)  to 
obtain  the  hydrodynamic  parameters.  To  minimize  the  error,  define  the  scalar  positive 

squared  error  measure,  J(n)  =  f,e'{t)e(t)  then  the  minimization  ofJis  given  by 

i=\ 

^  =  0  =  -lH'(t)e{t) 

de 

and  substituting  from  Equation  (29)  yields 

(32) 

/=1 

Rearranging  Equation  (30)  to  get  the  generalized  equation 

iH'(t)y(t)=iH'{t)H(t^(t)  (33) 

t=l  /-I 

Writing  the  equations  in  matrix  form  and  solving  for  0(f)  gives 

0  =  [H'H^H'y  (34) 

Using  the  matrix  divide  function  in  MATLAB  to  evaluate  Equation  (32)  will  give 
a  result  that  is  a  least  squares  fit  for  0 (?) . 
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The  same  result  may  be  found  using  guassian  random  assumptions  in  which  the 
solution,  0  in  Equation  33  is  the  most  likely  solution  in  which  its  probability  is  maximum 

Gelb  (1974).  It  should  be  noted  that  the  regression  matrix,  [  Xi/'(t)tf(f)]  must  be  positive 

/=] 

and  strong  (with  no  singularity)  otherwise  its  inverse  does  not  exist.  This  means  that  the 
system  must  be  perpetually  excited  by  its  input. 
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IV.  DATA  ACQUISITION  AND  PROCESSING 


This  chapter  discusses  the  data  collected  by  the  ARIES  AUV  including  the 
sensors.  How  the  data  is  processed  and  the  procedure  used  to  obtain  the  hydrodynamic 
coefficients. 

A.  AUV  DATA  ACQUISITION 


An  onboard  computer  samples  all  the  data  obtained  from  the  AUVs  sensors  at 
eight  hertz.  This  data  is  maintained  onboard  until  the  AUV  surfaces  and  transmits  the 
data  through  an  aerial  link,  or  it  is  downloaded  by  use  of  an  interface  cable  connected  to 
the  AUV  after  recovery.  The  parameters  logged  by  the  AUV  are  post  processed  using  the 
MATLAB  script  file  in  Appendix  2.  All  the  parameters  logged  are  listed  below  by  sensor 
and  variable  name  assigned  to  it.  Some  key  variables  recorded  at  8  Hz.  are 

1.  Pressure  Cell 


a.  Depth 

2.  Gyro/Compass/Acoustic  Doppler  Log 

a.  Heading 

b.  Forward  Velocity  (Relative  to  ground  as  well  as  the  water  mass) 

c.  Sideslip  Velocity  (Relative  to  ground  as  well  as  the  water  mass) 

d.  Heave  Velocity  (Relative  to  ground  as  well  as  the  water  mass) 

e.  Roll  Angle 

f.  Roll  Rate 

g.  Pitch  Angle 

h.  Pitch  Rate 

i.  Yaw  Angle  (Heading)  from  Magnetic  Compass 

j.  Yaw  Rate 

k.  Altitude  above  bottom 

3.  Global  Positioning  Receiver  (GPS) 

a.  Position 

b.  Heading 
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c.  Speed 

d.  Differential  Correction  Status  (0/1) 

4.  Onboard  Computer  (Outputs  from  Controllers) 

a.  Rudder  Position  Ordered 

b.  Plane  Position  Ordered 

c.  Thruster  Command  Voltage 

B.  DATA  EXTRACTION  AND  PROCESSING 

The  raw  data  is  post  processed  using  MATLAB  file  assign38.m  shown  in 
Appendix  B.  The  assign38.m  file  also  produces  a  plot  of  vehicle  position  during  the 
mission.  Figure  5  shows  the  track  of  the  ARIES  AUV  mission  used  in  this  report.  For 
the  analysis  procedure  used  in  this  report,  a  set  of  continuously  changing  (dynamic)  data 
was  required.  This  data  set  represents  a  mission  to  track  a  set  of  waypoints  forming  a 
lawnmower  pattern.  On  leg  3  an  oscillating  response  produced  data  to  identify  the  AUV 
parameters. 
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Path  plot  X,Y 


Figure  5  MATLAB  Plot  of  Full  AUV  Mission  File  dl00600_04.d 


Figure  6  shows  the  subset  of  data  seen  in  Figure  5  that  is  used  in  the 
determination  hydrodynamic  parameters.  The  limit  cycling  behavior  arose  from  the  use 
of  an  initial  steering  controller-not  well  behaved-but  good  for  the  purpose  of  steering 
parameter  identification.  Once  this  portion  of  the  data  was  isolated  for  analysis  it  was 
analyzed  for  its  usability. 
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Path  plot  X,Y 


Figure  6  MATLAB  Plot  of  AUV  Mission  File  Portion  Used  In  Anal>  sis 

One  of  the  limitations  of  the  data  was  that  the  rudder  signal  was  the  command 
signal  and  not  the  actual  position.  The  ARIES  rudder  motors  are  stepper  motors  run  open 
loop.  Figure  7  shows  the  recorded  rudder  command  and  yaw  rate  for  the  portion  of"  the 
data  being  used.  Looking  closely  at  Figure  7,  it  can  be  seen  that  the  change  in  >  au  rate  in 
response  to  rudder  action  was  delayed.  This  delay  is  primarily  the  result  of  rudder 
response  time.  From  the  data  shown  in  Figure  7  and  observed  response  time,  it  was 
decided  to  use  a  first  order  lag  of  0.26  seconds  for  the  rudder  signal.  Figure  8  shows  one 
rudder  cycle  and  a  first  order  lag  version  of  the  command  signal.  A  first  order  lag  was 
used  to  more  accurately  describe  the  time  history  required  for  the  rudder  to  change 
position  in  response  to  the  command  signal.  To  reduce  bias  in  the  data,  all  parameters 
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were  converted  to  zero  mean  before  being  used  for  calculations.  The  data  is  now  ready  to 
be  used  to  calculate  the  hydrodynamic  parameters.  The  MATLAB  file  thesis_values.m  in 
Appendix  B  was  used  to  process  the  isolated  data  using  mass  matrix  values  calculated 
from  Chapter  II  and  Equation  (19)  developed  in  Chapter  D3. 


Rudder  Command  and  Yaw  Rate 


Figure  7  MATLAB  Plot  of  AUV  Rudder  Order  Command  and  Recorded  Yaw  Rate  Used 

For  Analysis 


29 


I _ i _ i _ i _ i _ i _ I 

40  45  50  55  60  65  70 

Time  in  Seconds 


Figure  8  MATLAB  Plot  of  One  Cycle  Of  AUV  Rudder  Order  Command  Original  Data 

And  Lagged  Version  Used  For  Analysis 
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V.  RESULTS 


A.  DISCUSSION  OF  RESULTS 


The  first  step  in  calculating  was  to  put  all  of  the  variables  into  a  set  of  consistent 
units.  For  this  portion  of  the  report  all  calculated  and  measured  values  were  either 
verified  or  converted  as  required. 

Then  using  Equation  19  and  the  assumption  that  there  is  no  cross  coupling  and 
therefore  the  off  diagonal  terms  of  the  mass  matrix  are  zero,  two  independent  differential 
equations  can  be  written  as  follows: 

=YVrVr+(Yr-mUoy  +  Ys8r  (35) 

{lZ2-N,y  =  N,vr+Nrr+Ns8r  (36) 

where  the  notation  as  a  function  of  time  has  been  omitted  for  clarity.  Rewriting  Equation 
36  explicitly  in  terms  of  time  and  using  IN  =(izz-Ni.)>  a=NVr,  JS  =  Nr,  and  y  =  Ns 

results  in 

=  «v,+  y$t  (37) 


Solving  Equation  37  for  rl+l 


yields  the  result 


_  ccA t  f  PA/  ]  yA /  _ 

'+1  IN  '  [IN  /  IN  ' 


(38) 


Now  using  the  theory  developed  in  Chapter  IH  and  relating  Equation  38  to 
Equation  29  and  3 1  then, 

9(0 = [9i(4  e2 (/),  03 (/)] ,  y(t)  =  rt+1 ,  and  Hit)  =  [v, ,  r, ,  5, ] 
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aAf 


where  0,  =  ,  02  =  — +1  and  03  =— .  Casting  these  results  in  the  form  Equation  33, 


yA  t 


IN 


IN 


IN 


§  = 


H’y 

H'H 


where  H  is  a  matrix  and  y  is  the  vector  r,+] .  Using  MATLAB  to  perform  the 


matrix  operations  gives  a  least  squares  fit  for  0  as  discussed  in  Chapter  III.  From  the 
values  calculated  for  0,,  02  and  03 ,  the  values  of  a,  (i  and  y  (NVr  ,  Nr  and  Ns )  are  readily 

determined.  Similar  analysis  was  done  for  Equation  34  to  obtain  the  sway  parameters. 


The  MATLAB  file  thesis_values.m  calculated  the  hydrodynamic  coefficients 
Yy  ,Yr,Y6,N6,Nv  and Nr.  The  values  obtained  for  these  coefficients  are  summarized  in 

Table  7. 


Table  7  Calculated  Hydrodynamic  Steering  Coefficients  of  the  ARIES  AUV 


Yr 

Y* 

K 

68.16 

406.30 

69.9 

-10.89 

-88.34 

-35.47 

Once  these  coefficients  had  been  calculated  they  needed  to  be  evaluated  for 
accuracy.  The  coefficients  and  model  developed  were  combined  and  a  simulation  was 
conducted  to  compare  with  the  data  obtained.  Figure  9  shows  the  original  yaw  rate  plot 
with  the  plot  predicted  by  the  model.  Figure  10  shows  the  original  sway  plot  with  the 
plot  predicted  by  the  model.  It  can  be  seen  that  in  both  plots  the  two  are  very  close  to 
being  the  same,  thereby  showing  the  calculated  coefficients  are  accurate. 
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Sideslip 


Figure  10  MATLAB  Plot  of  Measured  And  Predicted  Side  Slip 

B.  ERROR  MEASURES 

The  scalar  positive  squared  error  measure  J(n)  (discussed  in  Chapter  III  C)  for 
the  coefficients  obtained  was  6.2055  x  10'5  for  sway  and  5.3467  x  10'6  for  yaw  rate. 
While  this  is  a  small  number,  the  question  remains,  how  does  this  compare  to  the  data  in 
a  magnitude  sense.  To  determine  the  relative  error,  the  mean  square  error 

N  i 

was  divided  by  the  square  of  the  data’s  standard  deviation  ay2.  Taking  the  square  root  of 
this  to  get  the  relative  error, 
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rel  = 


V 


N  t=\ 


xlOO% 


results  in  a  relative  error  of  4.38  and  1.79  percent  in  sway  and  yaw  rate  respectively. 
Thus,  the  maximum  error  overall  was  less  than  5  percent. 
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VL  CONCLUSIONS 


The  overall  problem  of  underwater  vehicle  modeling  is  very  complex.  To  obtain 
an  accurate  hydrodynamic  model  for  any  vehicle  that  is  not  a  standard  shape  by  purely 
analytical  means  is  extremely  difficult  at  best.  The  use  of  other  methods  such  as  use  of 
model  testing  in  tow  tanks  requires  additional  time,  travel  and  expense  and  do  not 
guarantee  accurate  results.  Thus,  the  method  used  in  this  report  has  many  positive 
features.  The  utility  of  this  approach  is  to  improve  autopilot  control  design  by  using 
engineering  approximations  to  the  model  of  a  vehicle.  The  data  collected  by  the  AUV 
can  then  be  used  to  determine  realistic  estimates  of  the  hydrodynamic  parameters.  An 
increased  knowledge  of  hydrodynamic  parameters  then  leads  to  a  better  control  system 
design  and  performance. 

While  this  approach  requires  that  an  initial  steering  controller  be  designed 
successfully  so  that  the  input-output  data  may  be  acquired,  the  results  of  this  study  has 
led  to  the  design  of  an  improved  controller  for  which  the  limit  cycling  shown  here  is  not 
present. 
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APPENDIX  A.  MAPLE  FILES  FOR  ARIES  MASS  MOMENT  OF 

INERTIA  CALCULATIONS 

This  appendix  contains  MAPLE  files  for  the  calculation  of  the  mass  moment  of 


inertias  for  the  ARIES  AUV  nose,  mid  and  tail  sections  respectively. 
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nose.mws 


%  Worksheet  to  calculate  AUV  nose  section  mass  moment  of  inertia. 
Calculate  Ixx 


>  restart; 

>  L :=37 . 5+Y : 

>  A : =2  7 : 

>  B :  =8  : 

>  C : =5 : 

>  dm:=pi*rho*A*C* (1-YA2/BA2)  : 
>x:=A* (1-YA2/BA2) A0.5: 

>  z :=C* (1-YA2/BA2) A0.5: 


Nose  section  equations 

Length  from  center  of  midsection  to  tip  of  nose  section 
Length  along  x-axis 
Length  along  y-axis 

Length  along  z-axis 
Differential  mass 
Equation  of  ellipsoid 


> m : =int (dm, Y=0 . . B  ) : 

> IAA:=pi/4*z*xA3*dy: 

>  dIAAm : =rho * I AA : 

>  dIx:=dIAAm/dy+YA2*dm: 

>  dlxa : =dIx*M/m : 


Calculate  mass  of  section 
IAA'  for  area  slice 
Thickness  =  dy 
Calculate  lx 

M  =  m  so  answer  is  in  form  M*Ixx 


> lx : =int (dlxa , Y=0 . . B) : 

>  simplify  (lx)  ;  Ixx 

158.6000000  M 


Calculate  Iyy 


>  restart; 

>A:=27: 

>B:=8: 

>  L : =37 . 5+Y : 

>  C : =5 : 

>  dm:=pi*rho*A*C* (1-YA2/BA2) 
>x:=A* (1-YA2/BA2) A0.5: 

>  z:=C* (1-YA2/BA2) A0.5: 

> m : =int (dm, Y=0 . .B  ): 

> IAA:=pi/4*z*xA3*dy : 

>  dIAAm : =rho *  IAA : 

>  dl Y :  =dIAAm/dy+LA2  *dm : 

>  dl Ya : =dl Y*M/ m : 

> IY:=int(dIYa,Y=0. .B) : 

>  simplify (IY) ; 


Differential  mass 
Equation  of  ellipsoid 

Mass 

LAA'  area  for  area  slice 
Thickness  =  dy 


Iyy 


1789.850000  M 
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Calculate  Izz 


>  restart; 

>  A :  =2  7 : 

>  L:=37 .5+x: 

>  B : =5 : 

>C:=8: 

>  dm:=pi*rho*B*C*  (l-xA2/AA2)  :  Differential  mass 

>  y :  =B*  (l-xA2/AA2)  A0 . 5 :  Equation  of  ellipsoid 

>  z :=C* (l-xA2/AA2) A0 . 5 : 

> m :  =int (dm, x=0 .  .A  ):  Mass 

>  IAA:=pi/4*z*yA3*dx:  IAA' area 

>  dIAAm :  =rho*  IAA :  Thickness  =  dx 

> dIZ : =dIAAm/dx+LA2*dm : 

>  dl Za : =dl Z  *M/m : 

>  int  (dIZa ,  x=0 .  .  A)  ;  Izz 

2316.425000  M 

> 
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midsection.mws 


%Worksheet  to  calculate  AUV  mid  section  mass  moment  of  inertia. 


>  restart; 
>A:=75: 
>B:=16: 
>C:=10: 


Length  along  x-axis 
Length  along  y-axis 
Length  along  z-axis 


> Ixx=l/12*M* (BA2+CA2) ; 

89 

Ixx  =  —  M 


> Iyy=l/12*M* (CA2+AA2) ; 


5725 


iyy= ~ylm 


> Izz=l/12*M* (AA2+BA2) ; 


M 


> 
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tail.mws 


%  Worksheet  used  to  calculate  AUV  tail  section  mass  moments  of  inertia. 


>  restart; 

> 

> 

> 

>  a :  =2  6 : 

>b:=16 
>c:=10: 

>  L :=x+37 . 5 : 

of  tail  section 

> dm:=rho*c* (-b/a*x+b) : 


>  m  :=  dm  dx 
•'o 

1  b 2  dm 


>  Iz  :=• 


12 


>  dlz  :=  Iz  +  L2  dm 


Mass*I 


>  Iz  :=  dlz  dx 

•'o 

> 

> simplify (Iz) : 

IzM 

m 

>  simplify  MIz) 


dx 


>  ly  :=  f  dly 

JyM 

m 

>  simplify(M/y) 

_  1  (b2  +  c2)  dm 

>&_:= — n — 

>  dlx  :=  lx  +  x2  dm 


Tail  Section  Equations 


Length  along  x-axis 
Length  along  y-axis 
Length  along  z-axis 

Length  from  center  of  midsection  to  face 

Elemental  mass 
Calculate  mass  of  section 

Calculate  Iz 


use  M=m  to  get  result  in  terms  of 


>  Iy_:=l/12*cA2*dm: 

>  dly  :=  Iy_  +  L2  dm 


Iz  result 

2190.250000  M 


Calculate  Iy 


Iy  result 
2177.250000  M 


Calculate  lx 
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c  :=  f  dlx 

•'A 


>  lx 

>  Mix  := 

> 

>  simplify(Mx) 


dx 

IxM 

m 


lx  result 
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APPENDIX  B.  MATLAB  FILES  FOR  DATA  ANALYSIS 

This  appendix  contains  the  MATLAB  script  file  for  analyzing  data  acquired  form 

the  ARIES  AUV.  The  first  file  assign38  converts  data  stored  in  the  AUVs  computer  to  a 
matrix  form  with  variable  handles.  The  second  file  auvcalc.m  determines  the  mass 
matrix  coefficients,  converts  units  from  English  to  SI,  checks  radius  of  gyration  and 
incorporates  modeling  approximations.  The  file  thesis_values.m  evaluates  a 
predetermined  portion  of  the  data  set  and  determines  the  hydrodynamic  parameters  based 
on  values  from  auvcalc.m. 
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assign38.m 


clear  all 
clc 

load  dl00600_04.d; 

d=dl00600_04; 

State  Id 


Nav_Id 

Month 

Day 

Year 

Hour 

Minute 

Second 

t 

X 

Y 

Depth 

Alt 

phi 

theta 

psi 

u 

v 

Depth_dot 

Alt_dot 

P 

q 

r 


d  ( : ,  1 )  ; 
d  ( : ,  2)  ; 
d  ( :  /  3 )  ; 
d  ( : ,  4 )  ; 
d(:,5)  ; 
d  ( : ,  6)  ; 
d(:,7)  ; 
d  ( : ,  8 )  ; 
d(:,9)  ; 
d  ( :  / 1 0 ) ; 
d ( : , 11) ; 
d ( : , 12) ; 
d ( : , 13) ; 
d ( : , 14  )  ; 
d ( : , 15) ; 
d  ( : ,  1 6 ) ; 
d  ( : ,  1 7 )  ; 
d(:,18)  ; 
d ( : , 19) ; 
d ( : , 20) ; 
d(:,21); 
d ( : , 22 )  ; 
d ( : , 23) ; 
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Bias_psi 

Bias_r 

Diff 

NSv 

LatDeg 

LongDeg 

TtTc 

n _ Is 

n_rs 

n_bvt 

n_sv  t 

n_blt 

n_slt 

delta_r 

delta_sp 


d  ( : ,  2  4 }  ; 
d  ( : ,  2  5 )  ; 
d  ( : ,  2  6 )  ; 
d  { : ,  2  7 )  ; 
d  ( : ,  2  8 )  ; 
d  ( : ,  2  9 )  ; 
d  ( :  /  3  0 )  ; 
d ( : , 31)  ; 
d( : , 32) ; 
d ( : , 33) ; 
d ( : , 34 ) ; 
d ( : , 35) ; 
d ( : , 36) ; 
d( : , 37)  ; 
d ( : , 38 ) ; 


dy=0 . 5 ; Yl= [ Y ( 1 ) :dy:491] ;  dx= (-19-X (1) ) /length ( (Yl)-rl) ; 

X1=[X(1) :dx:X ( 1) + (length (Yl) -1) *dx] ; 

Y2=[491:0. 5:691] ;X2=-19*ones (length (Y2) ) ; 

figure (1) , elf ,plot (Y,X, 'b. ' ,  Y1,X1, ' r . ’ , Y2,X2, ' r . ' ) , grid, zoom, 
equal; 

title ( '  Path  plot  X,Y,  blue.  Desired  Path  in  red'); 
figure (2) ,  plot (t, delta_r, t, r,  ' r . ’ ) ,  grid,  zoom 


axis 
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%analysis  of  model 


rhat=r; 


K=-0 . 1256; 

lamda=  0.2363; 

for  i=l : (length (t) -1) , 

rhat (i+1) =rhat (i) +0 . 125* (~lamda*rhat (i) +K*delta_r (i) ) ; 

end; 

figure (3) , elf ,plot (t, r,  'b. 1 , t, rhat, 'r.  ' , t , -delta_r . /2, 1 k. 1 ) , grid 


figure ( 4 ) , elf ,  plot (t, atan2 (v, u) ) 
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auvcalc.m 


%calculate  volume  of  AUV 

%  all  measurements  calculated  in  inches 

format  compact 

format  short  g 

clc 

total_wt^4 90/32. 17/12;  %  convert  wt  to  mass 

(lb  sA2  inA-l) 

sw_denisty=64/32. 17/12A4;  % (lb  sA2  inA-4) 

Msw_denisty=1025;  %kg/mA3 

batt_wt=23/32 . 17/12;  batt_num=6;  batt_vol=6*6*10; 

%  tail  section 
t_vert=10; 
t_lat-16; 
t_long=26; 

t_vol=l/2*t_vert*t_lat*t__long  %  2080  inA3 

%  Thruster  Vol  (negative) 

neg_vol=2*pi/4*10*5A2+2*pi/4*16*5A2;  %  1021.02  inA3 

neg_voly=2*pi/4*16*5A2; 

neg_volz=2*pi/4*10*5A2; 


%  center  section 

c_vert=10;  al=27.5;  a2=30.5;  bl=17.5;  b2=20.5; 

c_lat=16; 

c_long=75; 

c__vol__total=c_vert*c_lat*c__long 

c_vol=c_vert*c_lat*c_long-neg_vol 

c_voly=c_vert*c_lat*c_long-neg___voly 

c_volz=c_vert*c_lat*c__long-neg_volz 


%12000  inA3 
%1097 9  inA3 
%11372  inA3 
% 11607  inA3 


%  nose  section 
f illfactor=0 . 4 ; 

n__vert=10/2;  n_lat=16/2 ;  n_long-27/2; 

n__vol=.  5*4 /3*pi*njvert*n_lat*n_long* fill factor  %  452.4  inA3 


%  total 

full__vol=c_vol__total+n_vol/fillfactor+t_yol 

total_vol=t_vol+c_vol+n_vol 

calc_vol=  (total__wt+5/32 . 17/12)  /sw_denisty 

bouyancy 

vol_error= (total_vol-calc_vol) /total_vol*100 

homo_density=total_wt/total_vol 

inA-4 

density__w__o_battery=  (total_wt-batt_wt*batt_num) 
/(total  vol-batt  vol); 


%15211 

%  13511  inA3 
%  13365  51b  estimated 


%  9 . 4e-5  lb  sA2 


%moment  of  inertia  about  center  of  mid  section  UNIFORM 


c_Ix=89/3*homo_density*c_vol_total- . . . 

1/12  *homo_density* (3*2 . 5A2+10A2 ) - . . . 

l/12*homo_density* (3*2 . 5A2+16A2 )  %33 

c_Iy=5725/12*homo_density*c__vol_total .  .  . 
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-  ( l/12*homo_density* (3*2. 5A2 +10A2 ) +homo_density*bl A2) ...  %538 

-  (l/12*homo_density* (3*2 . 5A2+10A2) +homo_density*b2 A2 )  . . . 

-  ( l/2*homo_jdensity*2 . 5 A2+homo_density*alA2 ) . . . 

-  (l/2*homo_density*2 . 5 A2+homo_density*a2 A2 ) 

c_Iz=5881/12*homo_density*c_vol_total . . .  %552 

- ( l/2*homo_density*2 . 5A2+homo_density*bl A2 ) . . . 

- (l/2*homo_density*2 . 5A2+homo_density*b2 A2 ) . . . 

-  ( l/12*homo_density*  (3*2.5A2-i-16A2)  +homo_density*al A2 )  .  .  . 

-  (l/12*homo_density*  (3*2 . 5A2+16A2)  +homo_density*a2 A2 ) 


%moment  of  inertia  of  nose  about 
n_Ix=158 . 6*homo_density*n_vol 
n_Iy=1789. 85*homo_density*n_vol 
n_Iz=2316 . 4  25*homo__density*n_vol 

%moment  of  inertia  of  tail  about 

t _ lx—  427/3*homo_density*t_vol 

t_Iy=  2177 . 25*homo_density*t_vol 
t_Iz=  2190*homo_density*t_vol 


center  of  mid  section  UNIFORM 
%7 
%1 6 
%98 

center  of  mid  section  UNIFORM 
%28 
%  425 
£427 


%total  moment  of  inertia  about  center  of  mid 
Ix_total=t_Ix+c_Ix+n_Ix 
Iy_t°tal=t_Iy+c__Iy+n_Iy 
Iz  total=t  Iz+c  Iz+n  Iz 


section  UNIFORM 
% 68  lb  sA2  inA-l 
%1039 
£107  8 


%check  values  of  Ixx,  Iyy  &  Izz 

armx=sqrt (Ix_total/total_wt)  17 . 32  in 

army^sqrt ( Iy_total/total_wt )  %28 . 62 

armz==sqrt  ( Iz_total/total_wt )  12  9.15 


Iconvert  quantities  to  metric  units 


Mtotal_wt=4  90*0.4536 
Mlxx-lx__total*0 .0254*. 4536*9. 81 
MIyy=Iy_total*  0 .0254*. 4536*9.81 
Mlzz=lz_total*0. 0254* . 4536*9. 81 
Mvol_total=full__vol* .  0254  A3 
ckx=sqrt (MIxx/Mtotal_wt ) *100/2.54 
cky=sqrt (MIyy/Mtotal_wt ) *100/2.54 
ckz=sqrt (MI zz/Mtotal_wt ) *100/2 . 54 


%222.26  kg  wt  of  auv 
% 7 . 68  4  9  kg"mA2 
%  117 . 4 
%121 . 9 

%0. 24926  nA3 

%check  of  conversion  7.32 


%28 . 62 
%29 . 1 6 


%added  mass  calculations 

L=  64;  W=8;  H=5;  % inches 

%ecen-l-  (  (  (W+H)/2)/L  ) A2  £0.98  ave  dia 

ecen=l- ( (  7.53/L  )A2)  %const  vol  const  L  0.98616 

alpha=2* ( l-ecenA2) /ecenA3* (1/2* log (  (1+ecen) / (1-ecen) ) -ecen)  %0 . 08583 

beta=l/ecenA2-  (l-ecenA2) /2/ecenA3*log (  (1+ecen) / (1-ecen) )  %0 . 95709 

kx=alpha/ (2-alpha)  %0. 044839 

ky=beta/ (2-beta)  %0 . 9177 

kr=ecenA4* (beta-alpha) / (  (2-ecenA2) * (2*ecenA2- (2-ecenA2) * (beta-alpha) ) ) 
%0. 76391 


B_div_g=Mvol_total*Msw__denisty 

255.5 

Xudot=-kx*B_div_g 

Yvdot=-ky*B_div_g 

Nrdot=-kr*MIzz 


%kg  mass  of  water  vol  moved 

%kg  -11.456 
%kg  -234.47 
%kg  mA2  -93.13 
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Zwdot=Yvdot ; 
Mqdot=Nrdot ; 

Mma s  s  =Mt o t  a l_wt 
M Y=Mma s s - Y vdo t 
IN=MIzz-Nrdot 


%222.26  kg  wt  of  auv 
%456 . 73 
%215 . 04 
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thesis  values.m 


clear  all 
clc 

format  compact 
load  osciilation.mat ; 
%  d-os  dilation; 


State_Id 

= 

d(:, 

l) ; 

Nav  Id 

= 

d  ( : 

,  2) ;  Month 

= 

d  ( : ,  3 )  ; 

Day 

= 

d(:. 

4)  ; 

Year 

= 

d  ( : 

, 5 ) ;  Hour 

= 

d  ( : ,  6 )  ; 

Minute 

= 

d  ( : , 

7) ; 

Second 

= 

d  ( : 

,  8 )  ;  t 

= 

d  ( : ,  9)  ; 

X 

= 

d( : , 

10) ; 

Y 

- 

d  ( : 

, 11) ;  Depth 

= 

d ( : ,  12 )  ; 

Alt 

= 

d  ( : , 

13); 

phi 

= 

d  ( : 

,14);  theta 

= 

d ( : ,  15) ; 

psi 

= 

d(:, 

16)  ; 

u 

= 

d  ( : 

,17);  v 

= 

d ( : , 18 )  ; 

Depth__dot 

d  ( : , 

19); 

Alt_dot 

= 

d  ( : 

,20);  p 

= 

d  ( : ,  2 1 )  ; 

q 

= 

d  ( :  / 

22)  ; 

r 

= 

d  { : 

,23);  Bias  psi 

- 

d(:,24); 

Bias_r 

d  ( : , 

25); 

Diff 

= 

d  ( : 

,26);  NSv 

= 

d  ( : , 27 ) ; 

LatDeg 

d  ( : , 

28); 

LongDeg 

= 

d  ( : 

,29);  TtTc 

= 

d  { :  ,  3  0 )  ; 

n_ls 

d( : , 

31)  ; 

n  rs 

= 

d  ( : 

, 32) ;  n_bvt 

= 

d  ( : ,  33)  ; 

n__svt 

= 

d  ( : , 

34)  ; 

n_blt 

= 

d  ( : 

, 35) ;  n  sit 

= 

d  ( : ,  36)  ; 

delta_r 

= 

d  ( : , 

37)  ; 

delta_sp 

SB 

d  ( : 

,38)  ; 

%dy=0 . 5; Yl 

.= 

[Y(l) 

:dy: 

491];  dx= ( 

-19-X ( 1 ) ) /length ( (Yl)+1) 

; 

%X1=[X(1) :dx:X (1)  +  (length (Yl) -1) *dx]  ; 

%Y2= [491: 0.5: 691] ;X2=-19*ones (length (Y2) ) ; 

% figure (1)  , cl f .plot (Y.X) %,  'b. Y1,X1,  'r. ' ,Y2,X2,  Tr. 1 ), grid, zoom, 
equal; 

%title(*  Path  plot  X,Y,  blue.  Desired  Path  in  redT ) ; 

%  figure (2) .plot (Y (2000 : 5000, 1) , X (2000 : 5000 , 1 ) ) 

%figure (2) ,  plot (t, delta_r, t, r, ' r . f ) , grid, zoom 

^analysis  of  model 

%rhat=r; 

K=-4 ; 

lamda=  234 ; 

%for  i=l : (length (t) -1) , 

%  rhat ( i+1 ) -rhat ( i ) *lamda+K*delta_r ( i ) ; 

%end; 

%clf ,  plot (t, atan2 (v, u) ) 

MY=4  5  6 . 7  3 
IN=215 . 04 


rdr=delta_r*pi/180 ; 
rdr=rdr-mean (rdr) ; 

Rr=r*pi/180 ; 

Rr=Rr-*mean  (Rr )  ; 

Vv=v~mean  ( v)  ; 

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  lag  rudder 

y—  [  ] ; 


axis 
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y (1 ) =rdr (1) ; 
n=max (size (rdr) ) ; 

p=0 .2  6; 

q=0 .26; 
for  i=2:n 

y (i)  =  ( 1-p) *y (i-1 ) +q*rdr ( i-1 )  ; 

end 

yrdr=y 1 ; 

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 

%%%adjust  v  from  nose  to  center 

Lnose2center=55*2 . 54/100;  %mid  point  length  converted  to  meters 

Mid_Vv=Vv-Lnose2 center *Rr; 

%%%%%%%%%%%  cala  theat as %%%%%%%% % 

Hl=[Mid_Vv (1:1500, : ) , Rr ( 1 : 1500, :), yrdr ( 1 : 1500, :)]; 

Yl=[Mid_Vv(2:1501, : ) ] ; 
thetal=inv (HI T  *H1 ) *H1 ’ *Y1 

H2=[Rr (1:1500, : ), yrdr (1 : 1500,  :)]; 

Y2= [Rr (2 : 1501, :)] ; 
theta2=inv (H2 1 *H2 ) *H2 1  * Y2 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
rr=[]  ;  rr  ( 1 )  =Rr  ( 1 )  ; 

for  i=l : 1500 ,  rr (i+1 ) =rr (i ) +0 . 125* (-67 . 9/IN*rr (i) -32 . 4/IN*yrdr (i) ) ;end 

%%%%%%%%%%%%%calc  paramenters 
% 

Yv= ( thet al ( 1 ) -1 ) *  8  *MY 
Yr=thetal (2) *8*MY 
Yd=thetal (3) *8*MY 
%Nv=theta2 (1) *IN*8 
Nr=(theta2 (1)-1) *IN*8 
d=theta2 (2)*IN*8 

%check  MX 
det (HI 1 *H1 ) 
cond (HI ’ *H1 ) 

%M= [ IN] ; 

%A= [Yr] ; 

%MMM=MA - 1 *A 
%eig (MMM) 

%%%%%%%%%%%check  error 
errorl=Yl-Hl*thetal ; 

Jl-l/max (size (HI) ) * ( error 1 f  *errorl)  %6. 22e-5 

error2=Y2-H2  *theta2 ; 

J2=l/max (size (H2) ) * (error2  ’  *error2)  %5 . 41e-6 


%4507 
%422 . 36 


% [0.12268; 6. 9742] 


figure  (2 ) 

subplot (311) 

plot (t , Rr, t, yrdr,  1 : 1 ) 

xlabel(Tv  ...  corrected  for  position’) 
subplot (312) 
plot (HI ( : , 2 ) ) 
xlabel(’yaw  rate’) 
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subplot (313) 
plot (errorl, ' r . f ) 
xlabel ( 1  error  1 ) 

%figure (3) , plot (phi) 

%clf , plot (t, r,  '  b. 1 , t, rhat ,  '  r . ’ ,  t , -delta_r . 12,  ’ k. 1 ) f grid 
% figure (4 ) ,  plot (Depth) 


phidot=phi; 
for  i=2 : 1500 

phidot (i) = (phi (i+1) -phi (i-1) ) *4 ; 

end 


figure (6) 
subplot (312) 

plot (t,rr, 1 r . 1 ,t,Rr, ,bl ) 
subplot (311) 
plot  (errorl ) 
subplot (313) 

plot (r, 'r1 ) 
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